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Abstract. Gravitational wave observations can potentially measure properties of neutron star equations 
of state by measuring departures from the point-particle limit of the gravitational waveform produced in the 
late inspiral of a neutron star binary. Numerical simulations of inspiraling neutron star binaries computed 
for equations of state with varying stiffness are compared. As the stars approach their final plvmge and 
merger, the gravitational wave phase accumulates more rapidly if the neutron stars are more compact. 
This suggests that gravitational wave observations at frequencies around 1 kHz will be able to measure a 
compactness parameter and place stringent bounds on possible neutron star equations of state. Advanced 
laser interferometric gravitational wave observatories will be able to tune their frequency band to optimize 
sensitivity in the required frequency range to make sensitive measures of the late-inspiral phase of the 
coalescence. 



1. Introduction 

In numerical simulations of the late inspiral and merger of binary neutron-star systems, one commonly 
specifies an equation of state for the matter, perform a numerical simulation and extract the gravitational 
waveforms produced in the inspiral. The scope of this talk is to report work on the inverse problem: 
if gravitational radiation from an inspiraling binary system is detected, can one use it to infer the bulk 
properties of neutron star matter and, if so, with what accuracy? To answer this question, a number of 
simulations with systematic variation of the equation of state (EOS) is performed. The simulations use 
the evolution and initial-data codes of Shibata and UryQ to compute the last several orbits and the merger 
of neutron stars, with matter described by a parameterized equation of state, previously obtained in [1]. 
Our analysis uses waveforms from these simulations to examine the accuracy with which detectors with 
the sensitivity of Advanced LIGO can extract from inspiral waveforms an EOS parameter associated 
with the stiffness of the neutron star EOS above nuclear density. 

One might expect that gravitational wave observations can potentially measure properties of neutron 
star EOS by measuring departures from the point-particle limit of the gravitational waveform produced 
during the late inspiral of a binary system. A way to quantify this departure is to make use of the fact 
that the gravitational wave phase near the end of inspiral accumulates more rapidly for smaller values 
of the neutron star compactness (the ratio of the mass of a neutron star to its radius). In this way, 
greater accuracy may be obtained compared to previous work that suggested the use of an effective 
cutoff frequency to place constraints on the equation of state. Our results based on this approach suggest 
that realistic equations of state will lead to gravitational waveforms that are distinguishable from point 



particle inspirals at a distance < 100 Mpc for an optimally oriented system, using Advanced LIGO 
in a broadband or narrowband configuration. Waveform analysis, that uses the sensitivity curves of 
commissioned and proposed gravitational wave detectors, allows us to also estimate the accuracy with 
which neutron star radii, closely linked to the EOS parameter varied, can be extracted. The choice 
of EOS parameter varied in this work is motivated by the fact that, as Lattimer and Prakash observed 
[3], neutron-star radius is closely tied to the pressure at density not far above nuclear. Analysis of the 
numerical waveforms also indicates that optimizing the sensitivity of gravitational wave detectors to 
frequencies above 700 Hz can lead to improved constraints on the radius and EOS of neutron stars. 



2. Equation of state 

Difficulties inherent in astrophysically constraining the EOS of neutron star matter are the wide range 
of different microphysical models, implying the lack of a model-independent set of fundamental 
observables, and the fact that models often have more free parameters than the number of astrophysical 
observables. Therefore, an effort to systematize such constraints necessitates a form of EOS (i) 
parameterized in a model independent way, with a number of EOS parameters (ii) large enough to 
accurately model the universe of candidate EOS, but (iii) smaller than the number of neutron star 
observables that have been measured or will have been measured in the next several years. 

In the past, for different purposes, Miiller and Eriguchi [4] have used non-relativistic piecewise 
polytropes, with a large number (20-100) of segments, to accurately approximate realistic EOS, in order 
to construct Newtonian models of differentially rotating stars. More recently, with conditions (i), (ii) and 
especially (iii) in mind. Read et al [1, 2] have developed a relativistic piecewise poly tropic approximation, 
which reproduces the features of most realistic neutron star EOS with rms error < 4%, using only 3-4 
segments (or 3-4 free parameters), as described below. 



2.1. Astrophysical constraints on a piecewise polytropic equation of state 

The EOS pressure p is specified as a function of rest mass density p (rest mass density p = rribn is 
proportional to number density n, where m\, = 1.66 x 10~^^g is the mean baryon rest mass). The 
relativistic piecewise polytropic EOS proposed in [1, 2] has the form 

p{p) = Kip^^ (1) 

in a set of intervals < p < Pi'm rest mass density, with the constants Ki determined by requiring 
continuity on each dividing density pi. 

= (2) 
The energy density e as a function of p is then determined by the first law of thermodynamics, 

d- = —pd-, (3) 
P P 



which is integrated, to give 



6(p) = (l + a,)pc2 + ^ (4) 



(for Fj 7^ 1), where we used the condition e/p ^ as p ^ 0. This condition implies ao = and the 

other constants are fixed by continuity of e{p) at each pi. In the above equations, p is measured in 
g cm~^, e and p have units of erg cm~^ = dyn cm~^ and c is the speed of light. The specific enthalpy, [), 
is dimensionless and is given by 

f)=^ (5) 
pc^ 



The crust EOS is considered to be known and, for the purposes of this work, can be modelled with 
a single polytropic region, fitted to a tabulated crust EOS for the region above neutron drip (10^^- 
10^^ gcm~^); the numerical simulations considered do not resolve densities below this. This polytrope 
has Tq = Fcrust = 1.3569, with Kq = i^crust chosen so that p = 1.5689 x lO^^dyncm^^ when 
p = lO^^gcm^'^. The core EOS is constructed independently of crust behavior. In [1] it was found 
that three zones within the core, with adiabatic exponents Fi, Ts as parameters, are needed to model 
the full range of proposed EOS models with sufficient accuracy. A fourth parameter, pi = p{pi), is 
needed to fix an overall pressure shift, specified at the fiducial density pi . The dividing densities pi , p2 
are not considered new parameters. Instead, they are fixed by minimizing the error in fitting a large set 
of candidate EOS, leading [1] to the preferred values pi = 10^^'^gcm^'^ and p2 = lO^^ *^ gcm~^. The 
first dividing density po between the crust and core varies by EOS, but is also not a new parameter. It is 
defined as the density for which the log P vs. log p curves of the crust and core EOS intersect. Fixing 
the crust EOS and determining the three dividing densities in this way, allows one to model a large set of 
candidate EOS with only four parameters {pi, Fi, F2, F3}. 

Using this parameterized EOS, one can construct relativistic models of neutron stars and compute 
astrophysical observables, such as maximum mass, radius, moment of inertia etc. Doing so, while 
systematically varying the EOS parameters {pi, Fi, F2, F3}, Read et al [1] have obtained constraints 
on the allowed parameter space set by causality [Fig. 1] and by present and near-future astronomical 
observation. They find that two measured properties of a single star can potentially place stringent 
constraints on the EOS parameter space and that a potential measurement of the moment of inertia 
of a pulsar can strongly constrain the maximum neutron-star mass [Figs. 1, 2]. A detailed study of 
these constraints may be found in [1,2]. The potential combination of such astrophysical constraints 
with gravitational wave observation can restrict the allowed EOS parameter space region to a surface 
(thickened by the error bars of the observations) or even a line. The rest of this talk focuses on constraints 
placed on neutron star EOS by potential detection of binary inspiral waveforms. 




Figure 1. Left: Causality constraints are shown for a fixed value of Fi. For each EOS, the maximum 
speed of sound up to the central density of the maximum mass star is considered. The shaded surface 
separates the EOS parameter space into a region behind the surface allowed by causality, labeled Vg < c, 
and a region in which corresponding EOSs violate causality, labeled Vg > c. The second, outlined surface 
shows a weaker constraint to accommodate the expected error (~ 2 standard deviations a) in the speed 
of sound associated with a piecewise polytropic approximation to an EOS. 

Right: Surfaces representing the set of parameters resulting in a constant maximum mass. A single 
observation of a massive neutron star constrains the equation of state to lie above the corresponding 
surface. The exponent Fi is set to the least constraining value at each point (figures adopted from [1]). 



Figure 2. Left: Surfaces that represent the set of parameters that result in a star with a mass of 
1.338 M0 and a given moment of inertia / (from possible near-future observations of PSR J0737- 
3039A). Observation of moment of inertia / for a known mass M constrains the equation of state to 
lie on the corresponding surface. 

Right: Surfaces representing the set of parameters that result in a star with a mass of 1.4 Mq and a fixed 
radius R. Measurement of the radius i? of a neutron star with known mass M constrains the equation 
of state to lie on the corresponding surface. The wedge at the back right in both figures corresponds to 
incompatible values of pi and Fi (figures adopted from [1]). 

2.2. Initial candidate equations of state 

Systematic variation of the EOS parameters also allows us to determine which properties significantly 
affect the gravitational radiation produced, and thus can be constrained with detection of sufficiently 
strong gravitational waves. 

The models selected for this study use a variation of one EOS parameter in the neutron star core. We 
vary the core EOS by an overall pressure shift pi (specified at the fiducial density pi = 10^^'^ g cm~^), 
while holding the adiabatic index in the core regions fixed at Fi = r2 = Fa = 3. While only a 
subset of realistic EOS are well-approximated by a single core polytrope, reducing the EOS considered 
to this single -parameter family allows us to estimate parameter measurability with a reasonable number 
of simulations. 

After fixing the core adiabatic index, increasing the overall pressure scale pi produces a family of 
neutron stars with progressively increasing radius for a given mass; the choice of parameter pi is partly 
motivated by the observation of Lattimer and Prakash [3] that pressure at one to two times nuclear 
density is closely tied to neutron star radius, with R oc p\^^ . The radius is less sensitive to variation of 
the adiabatic index in the neutron star core, for reasonable adiabatic indices [1], as suggested by Fig. 2. 

Future work will incorporate additional variations of the EOS within the core. This could involve 
additional models of varying adiabatic index around a fixed pi (with Fi = F2 = F3 in all core layers), 
as well as multiple piecewise-polytrope zones (with Fi 7^ F2 7^ F3) within the core or EOS parameters 
yielding neutron stars of the same mass M and radius R but different internal structure. Such work 
would yield insight into the relative size and correlation of effects on the orbital evolution due to the 
stellar radius and internal structure. 

The first set of EOS were chosen to "bracket" the range of existing candidates, seen in Fig. 3. These 
models are HB with pi = 10^^'^'^ dyncm"^, a standard EOS; 2H with pi = lO"^^'^*^ dyncm"^, a stiff 
EOS; 2B, withpi = lO^^-^^ dyncm'^, a soft EOS. 

3. Numerical simulations 

For each parameterized EOS considered, the late inspiral and merger of a binary neutron star system is 
simulated. The gravitational mass of each neutron star in the binary is fixed to 1.35Mq, an average value 
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Figure 3. Initial choices of EOS for numerical evolution compared to the set of tabled EOS considered 
in [1]. Candidates are labelled in order of increasing softness: 2H, KB, 2B. Reproduced with permission 
from [56] ©(2009) by the American Physical Society. 



for pulsars observed in binary systems [32, 33]. The significance of tidal effects in this configuration is 
expected to fairly represent tidal effects over the relatively narrow range of masses and mass ratios in 
astrophysical binary neutron star systems. 

3.1. Initial data 

Initial data is generated by constructing quasi-equilibrium sequences of irrotational neutron stars 
in a binary system, following the methods of [12,21,34-37], using the Isenberg-Wilson-Mathews 
formulation [38,39]. Since tidal spin-up in neutron star binaries is estimated to be negligible [5,6], 
our initial data and quasi-equilibrium sequences are constructed assuming irrotational flow fields and 
neglecting the neutron star spin. The parameterized EOS of Eq. (1) is incorporated in the code to solve 
for initial data with a conformally flat spatial geometry coupled to the fluid equations of neutron star 
matter. The baryon number of each star is equal to that of a spherical isolated star with gravitational mass 
M = 1.35 M0. The initial data for the full numerical simulation is taken from the quasi-equilibrium 
configuration at a separation such that approximately 3 orbits remain until merger. Relevant quantities 
of the initial configurations for each parameterized EOS are presented in Tables 1, 2. 

3.2. Numerical evolution 

The Einstein equations are solved using the Baumgarte-Shapiro-Shibata-Nakamura formulation [40, 41]. 
The conformal factor of the spatial metric is evolved following [42] and the resulting black hole 
spacetime is handled using the moving puncture method [43, 44]. The lapse and shift are evolved using 
a dynamical gauge condition as in [42], while the Einstein equations coupled to the fluid equations are 
solved using the numerical scheme also described in [42, 56]. 

During inspiral, the fluid evolution is essentially free of shocks, whence the cold parameterized EOS 
specified in Sec. 2 is used in the simulations. During merger, when shocks are developed, we include a 
hot EOS component with a thermal effective adiabatic index Fcs, as described in [46]. Shock heating in 
the merger can increase the thermal energy up to ~ 20-30% of the total energy [17]. The merger and 



Table 1. Properties of initial candidate EOS. These range from the "softest" EOS at the top, which results 
in a prompt collapse to a black hole upon merger, to the "stiffest" at the bottom. Model HB is considered 
a typical EOS. The pressure pi, which is the pressure at density pi = 5 x lO^^gcm"^, determines 
the polytropic EOS for the neutron star core; all candidates have F = 3. Radius R and compactness 
GM/(?R are those of a single isolated spherical (TOV) star of gravitational mass M = 1.35 Mq, with 
radius measured in Schwarzschild coordinates. The maximum neutron star mass iWmax is given in the 
fifth column. 



Model 


logioPi [dyncm"^] 


ii[km] 


GM/c'R 


M^^ [Mq] 


2H 


34.90 


15.2 


0.13 


2.83 


HB 


34.40 


11.6 


0.17 


2.12 


2B 


34.10 


9.7 


0.21 


1.78 



Table 2. Quantities of initial data sets for irrotational binary neutron stars. Each star has a baryon 
number Mq equal to that of an isolated star with gravitational mass M = 1.35 M©. The ADM mass 
Madm of the initial slice includes the binding energy and J is the total angular momentum of the initial 
slice. 



Model 


Pmax [g Cm"^] 


Mo [Mq] 


MadM [Mq] 


cJ/{GMl^^) 


[Hz] 


2H 


3.73196 X 10^"^ 


1.45488 


2.67262 


0.993319 
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HB 


8.27673 X 10^4 


1.49273 


2.67290 


0.995361 


309.928 


2B 


1.38300 X 10^5 


1.52509 


2.67229 


0.987681 


321.170 



post-merger evolution is significantly affected by this hot component, since a collapse to a black hole 
may be delayed by the thermal pressure developed during and after merger. 

Gravitational radiation is extracted both by spatially decomposing the metric perturbation about flat 
spacetime in the wave-zone with spin-2 weighted spherical harmonics and by calculating the outgoing 
part of the Weyl scalar ^4. For equal mass neutron stars, such as those studied here, the quadrupole 

= 2,m = ±2) mode is much larger than any other mode, so we consider only this mode in this 
analysis. The waveforms output from the simulations are the cross and plus amplitudes h^c^D/GMioi 
and hxC^D/ GM^^i of the quadrupole waveform, as would be measured at large distance D 3> GMtot / 
along the z axis perpendicular to the plane of the orbit, versus the retarded time tret- Here Mtot is the sum 
of the two neutron star masses when they are far apart, Mtot = 1.7 Mq. The early part of the extracted 
waveforms, for t^st ^ 0, contains spurious radiation which is discarded in the data analysis. 

Each simulation is typically performed for three grid resolutions. In the best-resolution case, 
the diameter of each neutron star is covered with 60 grid points. Convergence tests with different 
grid resolutions indicate that, with the best grid resolution, the time duration in the inspiral phase is 
underestimated by about 1 orbit. This is primarily due to the fact that angular momentum is spuriously 
lost by numerical dissipation. Thus, the inspiral gravitational waves include a phase error, and as a result, 
the amplitude of the spectrum for the inspiral phase is slightly underestimated. However, it is found that 
the waveforms and resulting power spectrum for the late inspiral and merger phases, which we are most 
interested in for the present work, depend weakly on the grid resolution. 

4. Analysis of inspiral waveforms 

From the quadrupole waveform data, we construct the complex quantity 



h = h+ — ihx 



(6) 




Figure 4. Binding energy := Madm — -^tot as a function of orbital frequency fl for constant 
rest mass quasi-equilibrium sequences. The gravitational mass of eacii neutron star in isolation is 
M = Mtot/2 = 1.35 Af0. Departure from tlie point particle (3PN) curve occurs earlier for the 
stiffest EOS (candidate 2H - lowest compactness) and later for the softest EOS (candidate 2B - highest 
compactness). 



The amplitude and phase of this quantity define the instantaneous amplitude \ h\ and phase (j) = arg h of 
the waveform. The instantaneous frequency / of the quadrupole waveform is then estimated by 



The numerical data can be shifted in phase and time by adding a time shift r to the time series points and 
multiplying the complex h by e^^ to shift the overall wave phase by (/). 

It is useful to define a reference time marking the end of the inspiral and onset of merger. A natural 
choice for the end of the inspiral portion, considering the behaviour of a PP inspiral waveform, is the time 
of the peak in the waveform amplitude \ h\. However, the amplitude of the numerical waveforms oscillates 
over the course of an orbit. A moving average of the waveform amplitude over 0.5 ms segments is used 
to average this oscillation; the end of the inspiral is then defined as the time at the end of the maximum 
amplitude interval. The resulting merger time will be marked by solid vertical lines in the plots to 
follow. 

The numerical waveforms begin at different orbital frequencies. To align them for comparison, they 
are each matched in the early inspiral region to the same post-Newtonian point-particle (or PP) waveform. 



4.1. Post-Newtonian point particle 

Point-particle inspiral is not well-defined in full general relativity (GR), and one is left with the post- 
Newtonian point-particle (PP) approximation and fully general relativistic black-hole numerical solutions 
as natural substitutes. Fortunately, the Taylor T4 3.0/3.5 post-Newtonian specification, introduced in 
[47], agrees closely with numerical binary black hole waveforms for many cycles, up to and including 
the cycle before merger (see also [48]). This empirical agreement allows us to adopt the Taylor T4 
waveform as an appropriate PP baseline waveform, compatible with full GR until the last cycles. Quasi- 
equilibrium sequences [Fig. 4] and time-frequency analysis [Fig. 6] show that the binary neutron star 
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Figure 5. Solid lines show numerical waveforms, scaled by c^D/GMiou and aligned in time and phase 
to the same point-particle post-Newtonian inspiral (dashed line), using a method described in [56]. The 
two dashed vertical bars indicate the portion of the waveform used for matching; the last vertical bar 
indicates the end of inspiral time tu for the numerical waveform. The top two simulations, 2H and HB, 
show the start of post-merger oscillations from a hypermassive neutron star remnant in the simulation. 
2B shows quasinormal ringdown from a prompt collapse to a black hole following merger. Reproduced 
with permission from [56] ©(2009) by the American Physical Society. 

waveforms depart from this waveform 4-8 cycles (200-560 Mtot) before the best-fit PP merger, due to 
finite size effects. 

The TaylorT4 waveform is constructed by numerically integrating to obtain a gauge invariant post- 
Newtonian parameter related to the orbital frequency observed at infinity Q [56]. The orbital frequency 
evolution 0(i) and orbital phase evolution ^{t) are computed to 3.5 post-Newtonian order following 
[47]. For the numerical integration, one needs to specify the constants of integration by fixing 
coalescence time t^^ and the orbital phase at this time ^{t^^) = ^c^- These two parameters uniquely 
specify the 3PN waveform for given particle masses. The amplitude of the (l = 2,m = ±2) quadrupole 
waveform is then calculated to 3.0 post-Newtonian order as described in [49, 56]. 

To match the numerical data to the PP inspiral waveform, the two parameters, t^^ and are varied 
and the best match is obtained [56], to fix a relative time shift and a relative phase shift. The masses 
of the point particles in the PP waveform are fixed to be the same as the neutron stars in the numerical 
simulations (the gravitational mass of isolated spherical neutron stars with the same number of baryons) 
and so masses are not varied in finding the best match. With a goal of signal analysis, we choose the 
time and phase shift by maximizing a correlation-based match between two waveforms. More details on 
the matching technique may be found in [56]. We note, however, that longer simulations are required to 
more precisely fix a post-Newtonian match. 




4.2. Comparison of gravitational waveforms 

Unlike the case of matching binary black hole simulations to point particle post-Newtonian[47,48], 
the binary neutron star simulations show departure from point particle many cycles before the post- 
Newtonian merger time. Fig. 5 shows the four waveforms shifted so the best-match PP waveforms have 
the same and 0^^. As the stiffness of the EOS and thus the radius of the neutron stars, increases, the 
end of inspiral for the binary neutron stars is shifted away from the end of inspiral for post-Newtonian 
point particle. 
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Figure 6. Time-frequency behavior, vertical line markings as previous figure. The departure from the 
point particle time-frequency relations, shown using a long-dashed line, occurs between 700-1000 Hz 
depending on the EOS. In agreement with intuition and with Fig. 4, the departure occurs earlier for the 
stiffest and later for the softest EOS. Reproduced with permission from [56] ©(2009) by the American 
Physical Society. 

This can also be seen by plotting the instantaneous frequency of the numerical simulation waveform 
with the same time shifts, as in Fig. 6, which also shows more clearly the difference in the post-merger 
oscillation frequencies of the hypermassive remnants, when present. The larger neutron star produced 
by the stiff EOS 2H has a lower oscillation frequency than that from the medium EOS HB. The remnant 
forms with a bar-mode oscillation stable over a longer period than the ~ 10 ms simulated. The signal 
from such a bar mode may be even stronger when the full lifetime is included. Information that can 
be extracted from the presence (or absence) and characteristics of a post-merger oscillation signal from 
the hyper-massive remnant would complement the information present in the late inspiral. Although 
such information could probe the equation of state at densities higher than those present in ordinary 
neutron stars, aspects of the physics neglected in this study will likely come into play. Such aspects may 






include finite temperature, magnetic field, neutrino cooling and other effects that need to be accurately 
modeled, thus complicating the parameter extraction. As previously noted, our simulations included a 
hot component in the equation of state, which strongly affects the post merger behavior (and prompt 
versus delayed collapse to a black hole in particular), but inclusion of other effects is a subject for future 
study. 



4.3. Spectrum of gravitational waveforms 

Given h+ or , one can construct the discrete Fourier transforms (DFTs) /i+ or /ix . Both polarizations 
yield the same DPT amplitude spectrum \h\, with phase shifted by 7r/2, if one neglects discretization, 
windowing, and numerical effects (including eccentricity). The amplitude spectrum |^| is independent 
of phase and time shifts of the waveform. 

The stationary phase approximation is valid for the post-Newtonian waveform up to frequencies of 
about 1500 Hz (with < 10% error), so is used to plot the amplitude of the point particle spectrum. In 
terms of an instantaneous frequency 

the Fourier transform of the waveform has the amplitude 

For binaries comprised of equal mass companions, the gravitational radiation is dominated by quadrupole 
modes throughout the inspiral, so the analysis includes this mode only, as mentioned earlier. The wave 
phase (f) is negligibly different from twice the orbital phase 2$ until the onset of merger at very high 
frequencies so we can use the relation 

df '^d^ 
at TT at 

to write the amplitude of the Fourier transform entirely in terms of the functions dfl / dt and the amplitude 
^(/) = 1^1 of the polarization waveforms. 

The translation of emitted waveforms into the strain amplitude measured at a detector involves 
transformations incorporating the effects of the emitting binary's angle of inclination and sky location. 
These effects are absorbed into an effective distance -Deff, which is equal to the actual distance D for a 
binary with optimal orientation and sky location, and is greater than the actual distance for a system that 
is not optimally oriented or located. The detector will detect a single polarization of the waveform, some 
combination of the plus and cross polarizations of the emitted waveform. The polarizations extracted 
from the simulation can be used as two estimates of the strain at the detector for a given numerically 
modelled source, which give very close results and are subsequently averaged. 

To compare with noise curves in the usual units, the quantity f^^'^\h(f)\ is plotted, at a reference 
distance of Dgft = 100 Mpc and rescaling from previously plotted numerical output h{t)c?D /GM^oy 
using Mtot = 2.7 Mq. 

The full spectra of models 2H and HB, seen in Fig. 7, show peaks at post-merger oscillation 
frequencies. Waveforms of 2H and HB are truncated while the post-merger oscillation is ongoing; if 
the simulations were allowed to continue, these peaks would presumably grow further. The simulation 
of the softest candidate 2B, in contrast, promptly collapses to a black hole upon merger and has a short 
Uved, and relatively small-amphtude, quasinormal mode ringdown. 

Note that time-frequency plots like the ones in Fig. 6 showed that numerical waveforms follow the PP 
waveform at instantaneous frequencies of up to 700-1000 Hz, depending on the EOS. The disagreement 
in the spectra in Fig. 7 from the PP stationary-phase approximation waveform at frequencies below this 
is primarily due to the finite starting time of the numerical waveforms. To estimate spectra from the full 
inspiral, we construct hybrid waveforms. The short-term numerical waveforms are smoothly merged on 
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Figure 7. Colored curves show the DFT of full numerical inspiral waveforms, truncated and smoothly 
joined to best-match PP inspiral time series, at an effective distance Dgff = 100 Mpc. The stationary- 
phase point particle is shown by a dashed line for reference. Thick grey lines represent the noise spectra 
for Advanced LIGO (labelled "AdvLIGO" for the standard configuration and "Broadband" for the broad- 
band configuration) and the Einstein Telescope (labelled "ET"). Reproduced with permission from [56] 
(c)(2009) by the American Physical Society. 



to long-inspiral PP-PN waveforms, using Hann windowing to smoothly turn on or turn off a signal. To 
construct hybrid waveforms, these windows are used to turn on the numerical waveform as the matched 
post-Newtonian is turned off, such that the sum of the two window functions is 1 over the matching 
range. Because, as mentioned earlier, the dependence on the cold EOS is less straightforward during the 
merger, due to unmodelled physics that become important at that time, we focus instead on the signal 
from the waveforms during the inspiral region only, turning off the waveforms after the end of inspiral, 
tu- 

5. Distinguishability and detectability 

The question for neutron star astrophysics is whether these differences in the gravitational waveform will 
be measurable. We consider the possibility of detecting EOS effects with Advanced LIGO style detectors 
in varying configurations. 

We use several detector configurations commonly considered for Advanced LIGO include tunings 
optimized for 1.4 Mq NS-NS inspiral detection ("Standard"), for burst detection ("Broadband"), and 
for pulsars at 1150 Hz ("Narrowband"). The sensitivity is expressed in terms of the one-sided strain- 
equivalent amplitude spectral density Sh{f) (which has units of Hz^^/^) of the instrumental noise in 
Advanced LIGO. We also consider a provisional noise curve for the Einstein Telescope [51]. We 
consider only a single detector of each type, rather than a combination of detectors, for a preliminary 
estimate of detectability. 

Given two signals hi and /12, and a noise spectrum Sh{f), we define the usual inner product [52]: 

{hi\h2) = 4:Re J — — 



Table 3. g^fi in standard (NS-NS detection optimized) noise x (lOOMpc/Defr) 



Model PP 2B HB 2H 

~^P 0.32 0.55 0.69 

2B 0.48 0.63 

HB 0.60 

2H 



Table 4. g^f[ in broadband (burst-optimized) noise x (100Mpc/I?eff) 
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PP 


2B 


HB 


2H 


PP 





1.86 


2.67 


2.89 


2B 
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This inner product yields a natural metric on a space of waveforms with distance between waveforms 
weighted by the inverse of the noise. The detector output s is filtered against an expected waveform h 
using {s\h). Then 

^ ^ 7W) ^^^^ 

is the optimal statistic to detect a waveform of known form h in the signal s. If the detector output 
contains a particular waveform h that is exactly matched by the template used, then the expectation 
value of g is the expected signal-to-noise ratio or SNR of that waveform: 

g = y/{h\h). (13) 

Given a signal hi to be used as a template, one can ask whether a known departure from this signal 
can be measured. Assuming the modified waveform h2 is known, the expected SNR of /12 — h\ is 
similarly 

^diff = V {h2 - hi\h2 - hi). (14) 

We consider the two signals to be marginally distinguishable [53] if the difference between the 
waveforms hi and h2 has g^^i > 1. 

To analyze the measurability of finite size effects, we consider differences between hybrid inspiral- 
only waveforms matched to point particle PN waveforms in the early inspiral. After the numerical 
waveforms have been matched to the same post-Newtonian point-particle inspiral, the signals will be 
aligned in time and phase. We then then compare the resulting waveforms to each other, and to a 
PN-only waveform, using different Advanced LIGO noise spectra. We report results in terms of the 
SNR measured by a single detector at an effective distance Dgff = 100 Mpc. These results may be 
straightforwardly scaled to any distance in the wave zone, since h and g^fi are inversely proportional to 
the distance D. 

The difference between waveforms due to finite size effects is not detectable in the NS-NS detection 
optimized configuration of Advanced LIGO for ~ 100 Mpc effective distances. 

However, in both narrowband and broadband the differences can be significant, and waveforms are 
distinguishable from each other and from the PP waveform. Note that this impUes that, even if the details 
of a NS-NS late inspiral signal are not known, the difference between it and a point particle waveform 



Table 5. in narrowband 1150 Hz noise x (lOOMpc/Deff) 



Model PP 2B HB 2H 

~^P 091 3.69 2.12 

2B 2.92 1.45 

HB 2.25 
2H 



should be measurable. The quantity q^^{[ between the observed waveform and a best fit point particle 
waveform, limited to differences at high frequency, may be useful in itself to constrain possible EOS 
independently of waveform details. 



5.7. Parameter extraction 

Since the early part of the inspiral is well described by a point particle post-Newtonian approximation, 
we will assume that EOS effects on the waveform impact the late inspiral only. For simpUcity, we assume 
that orbital parameters, such as Mtot, mass ratio r], point particle post-Newtonian coalescence time t^^ , 
and phase shift cj)^, are determined from the observations of the earlier inspiral waveform, with sufficient 
accuracy that their measurement uncertainty will not affect the accuracy to which the late inspiral effects 
determine the EOS parameters. These measurements would be made by a broad-band instrument, in 
which the signal-to-noise ratio is expected to be high (~ 40 at 100 Mpc) and measurement accuracy is 
expected to be good [52]. Inaccuracies in these measurement could lead to biases in the measured EOS. 
This will be an important aspect to assess when high-quality binary neutron star simulations with various 
masses become abundant. 

With a one-parameter family of waveforms sampled, we can estimate the accuracy to which this 
parameter can be measured. There are also other EOS-related parameters which are not considered. In 
this first analysis we directly estimate the measurability of the EOS parameter pi , ignoring variations 
of r within the core. In the analysis of radius measurabihty, variations of the internal structures are 
respectively neglected. Expanding coverage of the EOS parameter space is underway. However, these 
initial parameter choices are expected to give the dominant contributions to finite size effects of the 
waveform. For example, the IPN tidal effect estimates of [7, 54] predict ~ 10% variation in the tidal terms 
contributing to binding energy and luminosity from changing internal structure — varying the apsidal 
constant — while keeping radius fixed. 

We estimate errors in parameter estimation to first order in 1 /~g or, equivalently, in S6^, using the 
Fisher matrix r^s = {dAh\dBh) [52]. Its inverse, (F"^)"^-^, yields 

§0A§0B = (r'Y^ (15) 
so that the expected error in a given parameter 9^ is 



{69^y = (r-i) 



UAA 



(16) 



and the cross terms of the inverse Fisher matrix yield correlations between different parameters. 

With a few simulations of varying parameter value, we estimate dh/dpi and dh/dR from two of the 
sampled waveforms, hi and /t2- For a single parameter 6 (which can be taken to be either pi or R), we 
have 



dh 



e={0i+e2)/2 



92 



hi 
9i 



(17) 



Table 6. Measurability of radius R (km) with broadband (burst-optimized) and narrowband (1150 Hz) 
configurations of Advanced LIGO, for a binary system at effective distance 100 Mpc. The values 
of R here are midpoint values of the radii of candidates shown in Table 1. The error 5R (km) is 
estimated using Eqs. (14), (18) (with g^ff taken from Tables 4, 5) and scales with effective distance 
as SRx (L»eff/100Mpc). 



Confij 


juration: 


Broadband 


Narrowband (1150 Hz) 


R = 


10.65 


±0.81 


±0.65 


R = 


13.40 


±1.52 


±1.60 



Table 7. Measurability of pressure pi (dyncm~^) at the fiducial density pi = lO^'^'^gcm^^, with 
broadband (burst-optimized) and narrowband (1 150 Hz) configurations of Advanced LIGO, for a binary 
system at effective distance 100 Mpc. The values of log^oPi here are midpoint values of logioPi from 
Table 1. The error (51ogioPi is estimated using Eqs. (14), (18) (with g^ff taken from Tables 4, 5) and 
scales with effective distance as Slogigpi x (Deff/lOOMpc). 





Model 


Broadband 


Narrowband (1150 Hz) 




rioPi = 34.25 


±0.13 


±0.10 




hoPi = 34.65 


±0.21 


±0.22 



where hi = h{9i) and ^2 = h{92), and then, for our one-parameter family where we neglect correlations 
with other parameters, we have to first order 



\2 



{h2 - hi\h2 - hi) 

Using adjacent pairs of models to estimate waveform dependence at an average parameter value, we 
then find estimates of radius measurability as shown in Table 6 and pi measurability as shown in Table 7 
for the burst-optimized noise and narrowband configurations. 



5.2. Sources of error 

Many higher order but likely relevant effects have been neglected in this prehminary analysis. Tidal 
effects may measurably influence the orbital evolution before the start of the numerical simulations, as 
estimated in [7], slowly enough not to be seen over the few cycles of the waveform matched to PP in 
this analysis. In one sense this analysis is a worst-case scenario, as it assumes exact PP behavior before 
the numerical match. Earlier drift away from point particle dynamics would give larger differences 
between waveforms, and more sensitive radius measurement, but poses a challenge by requiring accurate 
numerical simulation over many cycles to verify EOS effects. Combining numerical estimation with 
PN analyses hke those of [7, 8], and/or quasiequihbrium sequence information [Fig. 5] may clarify the 
transition between effectively PP and tidally influenced regimes. 

Some residual eccentricity from initial data and finite numerical resolution is present in the waveforms 
themselves. The eccentricity error may be estimated by comparing the plus polarization to the cross 
polarization shifted by 7r/2 from the same numerical waveform. This results in ^jj^g of ~ 0.3 for HB 
and 2H, rather than the expected quadrupole polarization cross-correlation of zero. The 2B waveform 
produces 'g^^ ~ 0.05 between polarizations. 

We have only a coarsely sampled family of waveforms; estimates of dhjdB are limited by this. 
Further, the length of the inspirals (see discussion on PN matching in [56]) Umits precision in choosing 



the best match time. We can estimate these effects on current results by varying the match region 
considered; this changes 'g^^ff by up to ~ 0.5 at 100 Mpc in the broadband detector. The resolution 
from existing numerical simulations is comparable to the difference between parameters of some of the 
models. The estimated uncertainty in each q^^^ is smaller than (but comparable to) the g^fi between 
candidates 2H, HB, 2B and PP. 

Finally, we note that use of a Fisher matrix estimate of parameter measurement accuracy is fully valid 
only in high SNR limit of Q^^ff > 10 [55], i.e., for distances < 20 Mpc in the broadband detector. The 
results do not take into account multiple detectors, nor multiple observations, nor parameter correlations. 
We conclude that, although these are of course first estimates, they should be better than order-of- 
magnitude. A full estimation of EOS parameter measurability will require more detailed analysis, with 
a larger set of inspiral simulations sampling a broader region of EOS parameter space, with mass ratios 
departing from unity, and with more orbits simulated before merger. 

6. Summary 

This is a first quantitative estimate of the measurability of the EOS of cold matter above nuclear density 
with gravitational wave detectors. Neutron stars of the same mass but different EOS and thus different 
radii will emit different gravitational waveforms during a late stage of binary inspiral. We calculated the 
signal strength of this difference in waveform using the sensitivity curves of commissioned and proposed 
gravitational wave detectors, and find that there is a measurably different signal from the inspirals of 
binary neutron stars with different EOS. 

Although the standard noise configuration of Advanced LIGO is not sensitive to differences in the 
waveform from finite size effects, a broadband (burst-optimized) configuration or a narrowband detector 
configuration of Advanced LIGO can distinguish neutron star EOS from each other and from point 
particle inspiral, at an effective distance < 100 Mpc. In addition, the provisional standard noise curve of 
the Einstein Telescope indicates the ability to resolve different EOS at roughly double that distance. 

With broadband Advanced LIGO at frequencies between 500 and 1000 Hz, our first estimates show 
that the radius can be determined to an accuracy of ~ 1 km x (lOOMpc/Deff)- Related first estimates 
show that such observations can constrain an EOS parameter parameter pi , the pressure at a rest mass 
density pi = 10^"^-^gcm~^, with an accuracy of 5pi ~ 5 x lO^^dyncm"^ x (lOOMpc/Deff)- These 
estimates neglect correlations between the parameters and other details of internal structure, which are 
expected to give relatively small tidal effect corrections. Although these results are preliminary, they 
strongly motivate further work on gravitational wave constraints from binary neutron star inspirals. The 
accuracy of the estimates will be improved with numerical simulation of more orbits in the late inspiral 
and a wider exploration of the EOS parameter space. 
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